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Abstract. 



Mode conversion in the solar atmosphere has been a problem of interest 



o . 

Recent numerical investigations of wave propagation near coronal magnetic null points 
, (McLaughlin and Hood, 2006) have indicated how a fast MHD wave partially converts into 

a slow MHD wave as the disturbance passes from a low-/? plasma to a high-/? plasma. This 
is a complex process and a clear understanding of the conversion mechanism requires the 
detailed investigation of a simpler model. An investigation of mode conversion in a strat- 
' ified, isothermal atmosphere, with a uniform vertical magnetic field is carried out, both 

numerically and analytically. Contrary to previous investigations of upward propagating 
waves (Zhugzhda and Dzhalilov, 1982a; Cally, 2001), this paper studies the downward 
propagation of waves from a low-/? to high-/? environment. A simple expression for the 
amplitude of the transmitted wave is compared with the numerical solution. 

o 

' Keywords: Sun: MHD oscillations, Sun: mode conversion 

1. Introduction 

> 

o 

\ for many years. It is important when the plasma /3, the ratio of the gas 

' pressure to the magnetic pressure, is close to unity; this is the region where 

. the sound and Alfven speeds are of comparable magnitude. Outside this 

\ area the modes are effectively decoupled, with one behaving like an acoustic 

■ wave and the other displaying a strong magnetic nature. Understanding the 

processes that occur as waves propagate through this /3 ~ 1 layer may in 
turn tell us the modes involved in coronal heating, and alert us to what we 
should be searching for observationally. 

Stein (1971) presented a full and detailed study of mode conversion across 
a density step, using the dispersion relations and boundary conditions to 
calculate transmission and reflection coefficients for fast, slow, and Alfven 
waves. The conditions under which the different modes are coupled to each 
other were also described. Zhugzhda (1979) and Zhugzhda and Dzhalilov 
(1981; 1982a,b) tackled the problem analytically, focusing on an isothermal 
atmosphere, stratified by gravity, and permeated by a uniform vertical mag- 
netic field. A solution was found in terms of Meier-G functions, specifically 
^24) allowing transformation coefficients to be found for waves propagating 
from high low in the atmosphere, to low /3, high in the atmosphere. This 
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study was extended by Cally (2001), who found that this solution could be 
expressed fully in terms of the equivalent hyper geometric 2-P3 functions. Re- 
flection, transmission and conversion coefficients were again found for waves 
incident from below, with another set of coefficients discovered in addition 
to those found by Zhugzhda and Dzhalilov (1982a). Our work uses the same 
simple one-dimensional set-up as these authors, but instead of studying the 
propagation from high /? to low (3 we examine propagation in the opposite 
direction. This study takes the form of a numerical investigation backed up 
by analytical work. Similarly, McLaughlin and Hood (2006) examined the 
behaviour of a fast magnetoacoustic wave travelling from a low-/? to a high- 
/? plasma, in this case as it approaches a magnetic null point. As the wave 
front crosses the (3 1 region, the fast wave can be seen to partially convert 
into a slow wave. 

Further numerical simulations have been carried out in two dimensions 
(Rosenthal et al., 2002; Bogdan et al, 2003), two-and-a-half dimensions 
(Khomenko and Collados, 2006), and in three dimensions (Carlsson and Bog- 
dan, 2006). All of these simulations found that the magnetic canopy (defined 
to be the isosurface where /3 = 1) is a crucial factor in mode conversion, along 
with the attack angle for the incident wave. This is supported analytically 
by Cally (2006), and by Bloomfield et al. (2006) who found observational 
evidence that the dominant wave mode changes from the fast mode to the 
slow mode as the waves propagate across the magnetic canopy. However, 
due to the complex nature of the numerical simulations it is difficult to 
disentangle individual physical effects. The simplicity of the model used in 
our investigation allows a much clearer picture to be built up, greatly aiding 
the understanding of the processes involved in mode conversion. This in 
turn will aid the interpretation of results from more complex models and 
simulations. 

The paper has the following outline. In Section 2 we introduce the MHD 
equations in the form that we shall be using them, detailing the assumed 
x-dependence and the non-dimensionalisation. The numerical simulations 
are discussed in Section 3, particularly in relation to their agreement with 
previous analytic studies. Finally we summarise our findings in Section 4. 



Throughout we work with the standard ideal MHD equations, namely 



2. Basic Equations 




V X (v X B) , 



- Vp xB)xB + pg, 
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|+V.(,v)=0, 

RpT 

P = 

where p is the mass density, v is the plasma velocity, p is the plasma pressure, 
fi is the magnetic permeability, B is the magnetic field, g is gravitational 
acceleration, 7 is the ratio of specific heats which we take to be 5/3, R is the 
universal gas constant, T is the temperature, and p. is the mean molecular 
weight which takes a value of 0.6 in the solar atmosphere. 

2.1. Basic Equilibrium 

A very simple, uniform, vertical equilibrium magnetic field, Bq = (0, 0, Bq), 
is chosen and the background pressure and density can be calculated from 
the momentum equation under equilibrium conditions of uniform tempera- 
ture (To) and the ideal gas law. 

Po{z) = Po(0)e-^/^, 
po (z) = po (0) e-'/", 

where H = RTq/ (Jig) is the coronal scale height 60Mm). Thus the 
plasma /? (= 2ppo/Bq) is dependent on z, ensuring that our waves will cross 
the (3^1 layer. 

To study the nature of this propagation we use the linearised MHD 
equations. Using subscript for equilibrium quantities and subscript 1 for 
perturbed quantities we have the linearised equation of motion 

Po^ = -Vpi + - (V X Bi) X Bo + Pig, (1) 
ot p 

the linearised induction equation 

= V X (vi X Bo) , (2) 
the linearised equation of mass continuity 

^ + V-(povi) = 0, (3) 
and the energy equation 

^ + (VI • V)po = — + (VI • V) poj . (4) 
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Henceforth it is assumed that we are working with the hnearised equations 
and subscripts on the perturbed quantities are dropped. 

2.2. x-Dependence and Non-Dimensionalisation 

To reduce our problem down to one dimension we assume that all variables 
are invariant in y, and choose the rc-dependence to be given by trigonometric 
functions of kx 

V = {vx (z, t) sin kx, 0, Vz {z, t) cos kx) , 

B = {Bx{z,t)smkx,0,Bz{z,t)coskx) , 

p = p {z,t) cos kx, 

p = p{z, t) cos kx. 

We also non-dimensionalise all variables by setting v = vqv, B = BqB, 
P = Poop, P = Poop, Po = PooPo, Po = PooPo, x = Lx, z = Lz, t = rt, and 
k = k/L, where a bar denotes a dimensionless quantity and vq, Bq, pqq, pqq, 
L, and r are constants with the dimensions of the variable they are scaling. 
We can then set Bq/^JJIpq = vq and vq = L/t, thus our speed is measured in 
units of Vq which represents a constant background Alfven speed. We also set 
(Sq = 2ij.poq/Bq = 2c|q/ (7^0) where Pq is the plasma /3 at a reference height 
of zero; note that g = c|q/ where c|q = 7Poo/poo- Thus we generate 

the dimensionless versions of Equations (1) - (4). Under these scalings t = 1 
(for example) refers to t = r = L/vq; i.e. the time taken for a wave to 
travel a distance L at the reference background Alfven speed. The bar on 
quantities is now dropped and it is understood that we are working with 
dimensionless values. 



2.3. Linearised Equations 

The linearised equations are 

1 dvx dBx kcl 

— -p + kBz 



v\ {z) dt dz 7 



1 dvz dp Lc^ 



0, 



vl{z) dt -1 dz i?7'^' 
dBx dvx 



dt dz 
dBz 
dt 



—kv. 



X, 



2 , .dp dvz L 

9^ + 57 = h""' - 

= 1^.-7/^^.. (5) 
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These may then be combined to form two second-order wave equations 

which are in agreement with Ferraro and Plumpton (1958). This form of 
the equations is much easier to analyse with the WKB approximation, as is 
done in Section 3.4. 



3. Numerical Simulations 

We numerically solve Equations (5) using the MacCormack method. This is 
a two-step, predictor-corrector, Lax-Wendroff method which is second-order 
accurate in time and space and, for linear harmonic waves, not strongly 
affected by numerical dispersion or diffusion. As we wish to study the prop- 
agation of waves incident from above, we drive Vz on the upper boundary, 
where we apply the conditions 

Vx = 0, 

Vz = siuLot, 



Bx 
dBz 

dt 
dp 

'd~z 
dp 

dz 



1 dB 



z 



k dz 

0, 
0, 
0. 



Imposing Vz on the upper boundary means that we are predominantly driv- 
ing a slow wave. Since the slow wave also has a small component of Vx, the 
condition Vx = means that there is a small component of the fast mode 
generated; however, this mode is evanescent and docs not propagate. To 
ensure that our code is as accurate as possible, we use backward differencing 
for the predictor steps and forward differencing for the corrector steps; thus 
we are using the corrected values on the upper boundary. The simulations 
are run for —10 < 2: < 5 and < t < 35; the end time is chosen so that 
we stop the simulation just before the wave front hits the lower boundary 
thus eliminating reflection effects. In all simulations, we choose Po = 0.2 
and L equal to the coronal scale height (H) so that z = 1 corresponds 
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to one coronal scale height. Having set these values, we are left with two 
free parameters: the driving frequency (a;) and the wavenumber in the x- 
dircction, k. By altering these parameters we can make comparisons between 
our results and different analytical models. 

3.1. Wave Properties 

As we wish to investigate wave propagation across the /? « 1 layer we 
shall describe the properties of waves above and below this region. An 
uncoupled slow magnctoacoustic wave (A; = limit) propagating through 
a low-/? plasma will change its behaviour to that of a fast magnctoacoustic 
wave as it passes into high-/? plasma. Similarly an uncoupled fast wave will 
change its behaviour to that of a slow wave as it travels from low to high-/? 
plasma. Despite this change in behaviour, the wave mode is the same, no 
conversion has occurred. Thus, when we discuss mode conversion, the slow 
wave driven on the upper boundary retains the properties of a slow wave 
as it propagates down into the high-/? plasma. We do not see any evidence 
of upward-propagating fast waves from the mode-conversion region. The 
transmitted part of the incident slow wave will continue into the high-/? 
plasma as a fast wave. 

This is demonstrated by the numerical simulations seen in Figure 1, which 
shows our variables under the transformations Vx u^e"^^^, Vz — >■ v^e^^'^, 
Bx — > Bxe'^/"^, Bz — s- -B^e"^/^, p pe~^/^, and p pe~^l'^. These transfor- 
mations remove the exponential behaviour, which is due to gravity, from the 
plots making other effects much clearer. The wave driven in on the upper 
boundary is, as discussed, predominantly a slow wave and this behaviour 
is retained until the wave reaches the Cg = va layer. At this point, mode 
conversion occurs, and in the high-/? plasma (to the left of the dashed line) 
the converted part of the wave propagates as a slow wave and the transmitted 
part as a fast wave. The plots of Vz,p and p in Figure 1 clearly show the fast 
wave out in front at z ~ —9.3, the interference seen for —4.9 < z < —1.8 is 
due to the interaction of the fast and slow waves. The plots for v^, Bx and 
Bz illustrate that the interference stops at z « —4.9 because this is as far as 
the slow wave has propagated. Gravity does not seem to affect the position 
of the mode conversion layer, which was also found to lie at Cg = in a 
weakly non-uniform medium. 

3.2. Varying the Wavenumber 

We now consider the effect of varying the wavenumber, k, on the mode 
conversion. Figure 2 shows the transformed vertical velocity, Vz, at t = 35 
for a range of values for k. The driving frequency, u = 2it, corresponds to 
driving in a wave with a period of one Alfven time, or 60s. 
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Figure 1. Transformed variables at t = 35, A; = tt, and uj = 27r. The Cs = va layer is 
denoted by a dashed line. 
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It is clear from the wave Equations (6) that when k = the fast and slow 
magnetoacoustic waves are completely decoupled. This can be observed in 
Figure 2 where no mode conversion is seen for k = 0, and indeed is not even 
visible for very small k. As k increases further the amount of mode conversion 
also increases, until almost none of the incident wave can be seen beyond the 
Cs = VA layer for A; = 7. In Section 3.3 we quantify the proportion of mode 
conversion using an approximation valid for small k, and in Section 3.4 we 
investigate behaviour for large k. 

3.3. Small k 

Cairns and Lashmore-Davies (1983, 1986) and Cairns and Fuchs (1989) 
developed a method of solving mode-conversion problems, which derives 
differential equations describing the coupled mode amplitudes from the local 
dispersion relation. These equations can then be solved analytically to find 
the transmission and conversion coefficients. In relation to our simulations 
this method can be applied when k is small, and to is sufficiently large 
compared to k. 

Starting with Equations (6) we take the time variation as e*'^* so that 
d/dt = iuj and let Vz = e^^'^Vz/cg. Under the assumption that k is small and 
uj large, the derivatives of and Vx with respect to z are larger than the 
functions themselves. Thus, our equations may be reduced to 

d iu) \ f d iu \ kCg dVz 

-r + — -; tt = 



dz VA {z) J \dz VA {z) J ^ VA {z) dz ' 

d f d iuj\ kcg dv^ 
+ — ]{- ]Vz 



^dz Cg J \dz Cs J VA (z) dz 

The equation for Vx is driven by the Vz that is imposed at the upper bound- 
ary. This inhomogeneous term has a wavenumber in z given by uj/cs and at 
zo, where Vz{zo) = Cg, there is a resonance and the amplitude of Vx increases 
rapidly while the amplitude of Vz is reduced. This is mode conversion. 

Expanding z = zq + around the mode conversion region {i.e. for the 
downward propagating waves described by the brackets containing minus 
signs) and letting d/dz = iu/cg away from this region, our equations may 
be written 

dvx . ( CO LJ A kc 



d^ \cs " 2va{z/'' 

dVz ii^^. kcs 
d^ cs 2vA {z) 



These satisfy the condition for energy conservation 



i(\v.\' + \Vzf)=0. 
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Eliminating Vx from Equations (7) and making the transformations given 
by Cairns and Lashmore-Davies 



we obtain 




(8) 



The solution of this equation is given by the parabolic cylinder function 
U{a,Q where 

ik'^Cs 1 
~ 2w 2' 

The asymptotic behaviour of these functions is described in detail in Abramowitz 
and Stegun (1964). 

Equations (7) may also be compared directly to the expressions used in 
Cairns and Lashmore-Davies (1983) to find the energy transmission and so 
the ratio of wave amplitudes in high and low-/3 is given by 

Thus, if we know the amplitude of the incident slow wave, we may multiply 
by A to predict the amplitude of the transmitted fast wave. 

A comparison between this analytical approximation and our numerical 
data is shown in Figure 3. The logarithm has been taken to highlight that 
the data is in excellent agreement for small k; no visible difference could be 
seen between the analytical prediction (9) and the original numerical data. 

Cally (2006) and Schunker and Cally (2006) apply ray tracing theory to 
the mode conversion problem. They also find these transmission coefficients. 

3.4. Large k 

For large k we may compare our numerical simulations with an analytical 
approximation derived by Roberts (2006) . In the limit of large k the first of 
our wave equations (6) reduces to 

(cj + v\ [z)) kvx + - ^) V, = 0. (10) 
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Figure 3. The logaxithm of the ratio of slow wave amplitudes in high and low (5. The solid 
line gives the values predicted using Equation (9) and the stars are values calculated from 
numerical simulations run for uj = A-n. 



We may use this to eliminate Vx from the first wave equation giving 



Vz = 0, 



where 1/cj, {z) = + l/v\ (z) is the tube speed. If we introduce 

\po (0)4(0); 

then our equation can be written in the form of the Klein-Gordon equation 



(11) 



^-4(^)g + ^^^(^)Q = o, 



(12) 



where 



n^{z) = c^{z 



2 / 1 4 (z) 4 (z) 



4 4 2c2?;2 [z] jv\ (z) \ c 



4(^) /^l 1 



7 



is a cut-off frequency. 
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Using the Klein-Gordon equation (12) with the term involving the cut-ofF 
frequency neglected, and the formal WKB method, as described in Bender 
and Orszag (1978), we find a solution of Q [z, t), valid for u large, to leading 
order 

r r I 



Q {z, t) ~ Dcj, (z) sin 



-dz 



Dc^ 



-2 



z } sm \ u) 
1 



i + — In 

Cs 

1 



ct{z) 

{CT (z) + Cs) {ct {Zrn) - Cs) 



{ct (z) - Cs) (ct (Zm) + Cs 



.Ct[Z) CT{Zm) 

where is a constant and Zm is the maximum value of z. We can then find 
an expression for Vz {z,t) using Equation (11) 

,1/2 



Vz {Z, t) 



VA jz) &/ {z„ 
VA {Zm) Or"^ (z) 
-2^ ' 



sm < LU 



i + — In 

Cs 



{cT {z) + Cs) {cT {Zm) - Cs) 



1 



(cT {z) - Cs) {cT {Zm) + Cs) 

(13) 

^Ct(Z) CT[Zm)J\) 

This expression agrees extremely well with the data obtained from the 
numerical simulations, as can be seen in Figure 4. An analytical approxima- 
tion to Vx can then be calculated using Equation (10); a comparison of this 




-10 -8 -6 -4 -2 2 4 

z/H 

Figure 4- Vertical velocity for uj — 4% and k = 407r at i = 35 (solid) and wave amplitude 
predicted from Equation (13) (dashed). 
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Figure 5. Horizontal velocity for uj = 4n and k — AOir at t = 35 (solid) and wave amplitude 
predicted using Equations (10) and (13) (dashed). 
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Figure 6. Transformed vertical velocity, Vz , for a; = A-k and k = 407r at t = 35 (solid) and 
wave amplitude predicted using Equation (13) (dashed). 
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result with the numerical data is shown in Figure 5. The agreement between 
the numerical simulations and the analytical approximation is clearly not as 
good for Vx (Figure 5) as it is for (Figure 4). If we remove the exponential 
behaviour from the vertical velocity, as is shown in Figure 6, we can see 
that for z < the agreement between the analytical and numerical results 
worsens. This is the same area for which we begin to see a deviation between 
the two solutions for the horizontal velocity. In order to improve agreement 
a higher-order WKB solution should be found. 



4. Conclusions 

In this paper we have considered the downward propagation of linear waves 

in an isothermal atmosphere permeated by a vertical background magnetic 
field. More specifically, we concentrated on the region where the wave passes 
from low-/? to high-/? plasma. As expected our simulations show mode con- 
version occurring as the wave passes through this region at the point where 
the sound and Alfven speeds are equal (Figure 1), similar to a weakly 
non-uniform medium. One of the free parameters in our simulations is the 
horizontal wavenumber (A;). Figure 2 shows that the value of k has a strong 
effect on the amount of mode conversion that occurs. For A; = 0, the two 
wave modes are completely decoupled and no mode conversion occurs, and 
as k increases so does the degree of mode conversion. As the wavenum- 
ber was shown to be an important factor, we then looked at analytical 
approximations in the limit of both small and large k. 

For small k we applied the method developed by Cairns and Lashmore- 
Davies (1983) for analysing mode conversion. This allowed us to find a 
solution for the vertical velocity [vz) which is valid in the mode-conversion 
region. This solution is related to the parabolic cylinder function, for which 
asymptotic solutions are readily available (Abramowitz and Stegun, 1964). 
We also found the amplitude of the transmitted wave, predicted from the 
amplitude of the incident wave, is in excellent agreement with our numerical 
data (Figure 3). This equation backs up our observation that the mode 
conversion increases with k and reveals that uj also has an effect on mode 
conversion, namely mode conversion decreases as lo increases. However, the 
effect due to the wavenumber k is the stronger of the two. 

For large k we followed an analysis carried out by Roberts (2006). From 
this we found a WKB solution for Vz valid for large uj (Equation (13)). This 
in turn was used to find an analytical solution for using Equation (10). 
Figure 4 shows very good agreement between the analytical and numerical 
results for the vertical velocity, however once the exponential behaviour is 
removed (Figure 6) we see that the analytical solution begins to deviate from 
the numerical results for z < 0. This is also true of the horizontal velocity 
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(vx) shown in Figure 5. Equation (10) linking to Vz itself is a very good 
approximation provided k is large. 

A thorough investigation of this very simple one-dimensional model has 
yielded some very interesting results, giving us some insight into the mode 
conversion problem. We hope to extend this work to include a non-isothermal 
atmosphere, making the model more physically realistic. 
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